epi draft 



Disordered Topological Insulators via (7*-Algebras 



o 

<N 

o 

o 

00 






^ 



O. 

o_ 

> 

m 

00 
00 

^, 

in 

o 
o 



Terry A. Loring^ and Matthew B. Hastings^ 

^ University of New Mexico - Department of Mathematics and Statistics, Albuquerque, NM 87131, 
1^ Microsoft Research - Station Q, Elings Hall, University of California, Santa Barbara, CA 
Department of Physics, Duke University, Durham, NC, 27708, USA 



USA 
93106, USA 



PACS 73 . 43 . Cd - Quantum Hall effects: theory and modeling 
PACS 72 . 25 . -b - Spin polarized transport 
PACS 73 . 20 . Fz - Weak or Anderson localization 

Abstract. - The theory of almost commuting matrices can be used to quantify topological ob- 
structions to the existence of localized Wannier functions with time-reversal symmetry in systems 
with time-reversal symmetry and strong spin-orbit coupling. We present a numerical procedure 
that calculates a Z2 invariant using these techniques, and apply it to a model of HgTe. This nu- 
merical procedure allows us to access sizes signiflcantly larger than procedures based on studying 
twisted boundary conditions. Our numerical results indicate the existence of a metallic phase in 
the presence of scattering between up and down spin components, while there is a sharp transi- 
tion when the system decouples into two copies of the quantum Hall effect. In addition to the 
Z2 invariant calculation in the case when up and down components are coupled, we also present 
a simple method of evaluating the integer invariant in the quantum Hall case where they are 
decoupled. 
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The study of topological insulators is one of the most ac- 
tive areas of physics today. Experimental and theoretical 
work has shown physical realizations of time-reversal in- 
variant insulators with strong spin-orbit coupling in both 
t^vo [1] and three dimensions [2 and a complete classi- 
fication of different insulating phases has been recently 
obtained using methods of Anderson localization [3] and, 
more generally, i^-theoretic techniques [3]. 

However, numerically it is difficult to determine the 
Z2 invariants that are signatures of topological insulating 
phases. For systems with translational invar iance, one can 
study the bundle over the momentum torus 0, while for 
systems without translation invar iance, Essin and Moore 
|6] were able to study the phase diagram of a graphene 
model by studying the model over a flux torus correspond- 
ing to twisted boundary conditions. Unfortunately, the 
flux torus approach is very computationally intensive: for 
each disorder realization, the Hamiltonian must be diago- 
nalized once for each point on a discrete grid on the flux 
torus, and then the connection on the torus must be com- 
puted. This limited the study to small systems, with at 
most 64 sites. 

In this paper, we present a different approach to calcu- 
lating a Z2 invariant, based on ideas in C* -algebras, in par- 



ticular the i^- theory of almost commuting matrices. We 
present a fast numerical algorithm based on these ideas. 
Computing the invariant requires a single diagonalization 
of the Hamiltonian, matrix function calculations on ma- 
trices at most half the size of the Hamiltonian, and finally 
the calculation of the Pfaffian of a real anti-symmetric ma- 
trix that is at most the size of the Hamiltonian. The most 
costly step is a single diagonalization, allowing us to study 
significantly larger samples, up to 1600 sites. 

Our invariant serves the same fundamental purpose as 
the invariant used in — proving that for certain low- 
energy bands it is impossible to find well- localized Wannier 
functions with time-reversal symmetry. These invariants 
are most likely equivalent, but that is another topic [7]- 

We apply this method to a model of HgTe [l], including 
the additional term i^BiA in [8, which breaks inversion 
symmetry. The Hamiltonian we use is ?^ of [8], with an 
additional on-site disorder term. This can be written as 
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i^BIA + V, 
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where in the notation of [8] 

/ -A\ 

^BIA = ^ (2) 

v-A ; 

is a /c-independent term that couples up and down spin 
components, and V is a random disorder term, diagonal 
in spin and band indices, described below. We map out 
part of the phase diagram, focusing on the transition with 
and without coupling between up and down spins. 

Band Projected Position Operators and Wan- 
nier Functions in Systems Without Time-Reversal 
Symmetry. — We describe first an integer invariant in 
the case of systems without time-reversal symmetry. We 
largely follow [9,, which described an integer invariant for 
two-dimensional systems on the sphere or torus, and a Z2 
invariant in the case of a sphere. In the next section, we 
present a Z2 invariant for systems on a torus. 

Consider a lattice Hamiltonian H for a single particle, 
tight-binding model on the surface of a torus. We assume 
that the hopping is short-range, so that the matrix element 
Hij is small if sites i and j are far separated on the torus. 
Let P be the projector onto the states below the Fermi 
energy Ep. If the Fermi energy lies in a spectral gap, 
then one can show using the locality of H that P is also 
local: Pij is small if i and j are far separated. In fact, P 
is still local even if the Fermi energy only lies in a mobility 
gap. We parametrize the position of a given site i on the 
surface of the torus by two angles, Oi and ^i, between and 
27r. Introduce two Hermitian matrices, O and ^. These 
are both diagonal matrices, with matrix elements B^^ = Oi 
and ^ii = (l)i. 

Using a basis of eigenstates of H, we can conju- 
gate the band-projected position matrices Pexp(z0)P and 
Pexp(i<l>)P by a single unitary to produce block matrices 



Pexp(ie)P 




U 



Pexp(z$)P 




V 



(3) 



where the first block corresponds to the empty states 
(those projected onto by 1 — P) and the second block cor- 
responds to the filled states (those projected onto by P). 

The matrices U and V are approximately unitary 
(U^U ~ / and V^V ~ /) whenever Ep lies in a mobil- 
ity gap and the system is large. To see this, note that if 
P is local, then P almost commutes with exp(iO). More 
precise results for a system on a sphere are in [9]. Further, 
U and V almost commute if Ep is in a mobility gap. 

Now we ask: given two matrices, U and F, which almost 
commute and which are almost unitary, does there exist a 
pair of matrices U\ V which exactly commute, are exactly 
unitary, and which are close to \J and Vl This question 
of approximating almost commuting matrices by exactly 
commuting matrices is an old problem in (7*-algebras [10]. 

Consider the quantity 
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tr(log(VC/y1'C/^')) = 27rim + r, 



(4) 



Fig. 1: Average index as a function of the Fermi level Ef-, for 
'w; = 20 (1200 samples), ^ = 30 (520 samples), ^ = 40 (800 
samples), it; = 60 (470 samples). Inset: scaling collapse with 
z/ = 7/3 [13]. 



where vn and r are real. We have exp(27ri7Ti + r) equal 
to det(FC/y^/7^). This determinant is real and positive 
(since det([/) = det([/t)) so vn is an integer. This integer 
m is precisely the topological invariant called the Bott in- 
dex. The approximation by exactly commuting matrices 
is possible if and only if tti = [12] . 

If U and V almost commute and are close to unitary, 
VUV^W is close to the identity. All of its eigenvalues 
stay away from the branch cut of the logarithm (chosen 
on the negative real axis), which is why the integer ?n is a 
topological invariant. If U and V exactly commute then 
m = 0. 

The physical importance of the Bott index is that when 
it is non-zero, it is not possible to find a complete, or- 
thonormal basis of localized functions (so-called Wannier 
functions) spanning the occupied states |9] . 

This invariant can be computed very quickly. We per- 
form a full diagonalization of the Hamiltonian to form the 
projector P onto states below Ef. We then construct U 
and V and calculate the trace of the logarithm in eq. (|1|) 
from the eigenvalues of VUV^W . The matrices U and V 
are of dimension equal to the number of occupied states. 
This is at most half the size of the dimension of H. 

Fig. [H shows the results of this invariant on a model of 
HgTe studied in [8j , with the up and down spin decoupled 
(A = 0). We chose constants A = 1, P = -1, P> = 0, M = 
—2, which sets the irrelevant terms of order k'^ to zero. The 
system had linear size w giving w'^ sites and 2w'^ states per 
spin component, and we added a diagonal disorder term 
on each site chosen uniformly from the interval [—4,4]. 
The transition sharpens as w increases. In the inset, we 
plot a scaling collapse with the exponent u = 7/3 from 
^13j . This collapse suggests that the transition is sharp in 
the thermodynamic limit. 

One can also use the commutator of U and y as a signa- 
ture of localization. In fig. O we plot this commutator as a 
function of Ep for various system sizes. One may see that 
as w increases, the commutator gets smaller away from 
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Fig. 2: Norm of commutator [f/, V] as a function of Ef- Inset: 
same plot, for A = 0.64, M = —1, with 300 samples for w = 40, 
and > 400 samples for smaller w. 



the transition. In the inset, we show the same plot for a 
system with A = 0.64, M = — 1, to show the existence of 
a metallic phase in this case. 

Finally, we mention the integer invariant on a sphere. 
To describe a sphere, we introduce three coordinate matri- 
ces Xi, X2, X3, labeling the positions of the sites in three- 
dimensional space, with the radius of the sphere set to 
one so that ^^ X^ = /. We then project these matrices 
into the band of occupied states to form three matrices 
which we call J^i, i^2, Hs. Then if Ep is in a mobility gap, 
^- Hf is close to the identity and Hi almost commutes 
with H^. We form the matrix 



BiHi,H2,H3) = J2H, 



(5) 



where a^ are the Pauli spin matrices. One may show that 
B^ is close to the identity if the Hi almost commute and 
^^ Hf is close to the identity. The integer invariant, also 
called the Bott index, is half the difference between the 
number of positive and negative eigenvalues of B. 

The reader may also ask about invariants on a system 
with open boundary conditions, such as a disk [14|, where 
there are two almost commuting band projected position 
matrices, PXP, PYP. In this case, there is no topological 
obstruction to approximating almost commuting matri- 
ces by exactly commuting matrices [T5l|l6] ; physically, we 
understand this as if the system on a disk is in a topologi- 
cally nontrivial phase then there will be gapless boundary 
modes and hence PXP, PYP will not almost commute. 

An Obstruction to Localized Wannier Functions 
with Time- Reversal Symmetry. — We now consider 
invariants in systems with time reversal symmetry. We 
begin with the invariant for matrices that almost repre- 
sent a sphere constructed in [9|. This theory applies to 
all triples {Hi, J^2, ^3) of self-dual Hermitian matrices for 
which 

\\[Hr,H,]\\<S, \\I-J2Hr\\<6 (6) 



with a small S. We need not convert a lattice over the 
torus to a lattice over the sphere; we compute the band- 
compressed periodic position operators and apply matrix 
functions to create matrices that satisfy eq. (|6|). 

Consider the matrix B{Hi, H2-, H^), which we define as 

1(1 + icjy (g) ay)B{Hi,H2,Hs){l - icfy x a^), (7) 

where Gy ® cFy is a product of the cFy matrices on the 
physical spin of the particle and on the pseudospin in- 
troduced to define B. The resulting matrix is real and 
anti-symmetric. The sign of its Pfaffian represents a po- 
tential Z2 obstruction in the group Z2 = { — 1,1}. 

We call this index the Pfaffian-Bott index. This index 
being —1 is first of all an obstruction to our being able to 
approximate the Hr simultaneously by exactly commuting 
self-dual matrices. The analysis in §V.B of [9] is valid in 
the self-dual case, and it shows that this is also an obstruc- 
tion to finding exponentially localized Wannier functions 
with time-reversal symmetry. 

The spectrum of B{Hi,H2^Hs) is pure imaginary and 
has a gap at zero. The smaller the commutators between 
the Hr the larger the gap, and the larger the distance must 
be to exactly commuting self-dual matrices. 

We prefer to study a physical system on a torus in order 
to have a regular lattice without defects. The Z2 invariant 
for the torus is not as simple as eq. (jl]). We compute the 
matrices U and V as above; these matrices /7, V will be 
self-dual if we are careful how we diagonalize the Hamil- 
tonian. We need to diagonalize via a symplectic matrix, 
found using a variation of an old algorithm J7] designed 
for Hamiltonians with time-reversal symmetry. 

We select a degree-one, continuous mapping from the 
torus to the sphere, a map (6^,0) ^ (^1,^2,^3) with 

XI = fW 

X2 = gW + Hcj)) cos 27ri0 (8) 

xs = h{(j)) sm27TiO 

We then apply this map to U and V by defining 

Hi = f{V) 

H2 = 9iV) + 1{HV),U^ + 1{HV),U} (9) 

Hs = l{h{V),U^}-l{hiV),U} 

and then compute the Z2 invariant of Hi,H2-, Hs . The 
anticommutators ensure that the H^ are self-dual. 

Expressing the torus invariant in terms of the 
sphere invariant requires computing with the matrix 
B{Hi, i^2, Hs)^ which has dimension twice as large as the 
number of occupied states, and so is at most the size of 
the Hamiltonian. We need only tridiagonalize the matrix 
B, which is real, further improving the speed. 

There are issues with applying the formulas in eq. ([9]). 
The matrix functions of V must be analytic functions, 
Laurent polynomials or polynomials in V and V^ since 
we don't expect V to be normal. That limitation is not 
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Fig. 3: Graphs of the /i, gi and hi as scalar functions 

compatible with defining a degree-one map from the torus 
to the sphere. This we overcome by choosing /, g and h 
that define a map to three-space that is close to the sphere. 
Our choice of / is similar that that used in |9l[Tl] , only 
with three derivatives vanishing at the extrema. See fig. [3l 
We define /i, gi and hi, truncating the Fourier series for 
gi and hi to produce order-5 trigonometric polynomials: 



/i(^) = ll sin(27r^) + ^ sin(67r0) + jfg sin(107r 




(10) 



Our invariant is most meaningful when ||1/''' — ]/ ^|| is 
small, so we replace gi{V) by a star-polynomial g{V) by 



9iiv) = e: 



bnV^' 



9iV) = E„=-5^"^'"^ where 



V^^^ equals V^ for positive n and {V^)~^ for negative n. 
This is much faster. Applied to scalars of unit modulus, 
there is no difference between gi{z) and g{z). We do a 
similar replacement of inverse by adjoint to define /(V) 
and h{V). The underlying map from the two sphere is 
illustrated by fig. E) 

Based on earlier experience [llj, we used also an alter- 
native method. We computed the polar decomposition 
F = Fo|V'| of F and a symplectic diagonalization WDW^ 
of Vb . This allowed discontinuous functions 



/((/.) = 1-20 
g{4>) = 

h{4>) = vi - f{4>r 



(11) 



as we can apply them to the diagonal of D. 

Both methods work well, with the star-polynomial 
method, eq. (fTQ|) . being faster. The log-based method, 
eq. (dU), gave smaller values for the important numbers 
||[i!f^,i!f5]|| and ||/ — X]^r||- The data in the plots were 
generated using the log-based method. These matrix func- 
tions, and 0{n^) algorithms for Pfaffians and symplectic 
diagonalization, will be discussed in [7j. 

The inability to approximate the H^ by self-dual, ex- 
actly commuting matrices implies the inability to approx- 




Fig. 4: Panel d) indicates using /i, gi and hi to map the torus 
(panel a)) close to the sphere. Panels b) and c) show the result 
of truncating the three functions to degrees one and three. 



imate the U and V by self-dual, exactly commuting ma- 
trices. We have an obstruction to the existence of expo- 
nentially localized Wannier functions with time-reversal 
symmetry in our original lattice model on the torus. 

Results for a System With Time-Reversal Sym- 
metry. — In fig. [5] we plot the average value of this Pfaf- 
fian as a function of Ef for a system using the model 
of HgTe studied in [8j. We include spin-orbit coupling, 
which in the terms of [8| means A = 0.64, M = —1. Un- 
like in fig. dj here the plots do not continue to sharpen as 
w increases, suggesting the presence of a metallic phase. 
When the U and V are far from commuting, the spec- 
tral gap of B{Hi,H2,Hs) is small and the value of the 
Pfaffian-Bott index is less meaningful. The wide plateau 
in fig. [2] inset is perhaps a stronger indication of a metal- 
lic phase. The fluctuating index implies that the phase is 
metallic. We are not certain of the interpretation of the 
index in the metallic regime, but we note that the index is 
a well-defined topological invariant so long as the Green's 
function has a sufficiently fast power-law decay. 

In fig. [5] inset, we compare for the case A = the average 
value of the Z2 invariant to the Z invariant for the same 
system, showing that the average value of the Pfaffian- 
Bott index is close to 1 — 2 * r?i. Recall that nonzero Bott 
index and Pfaffian-Bott index not 1 both indicate topo- 
logical insulators. 

In fig. [6l we show the phase diagram as a function of 
A with M = —1; in fact, A cannot be varied in experi- 
ment, but the ratio A/M can be as M depends on sample 
geometry [8]. For w = 20, we computed which values 
of Ef gave average index equal to ±1/3 or ±2/3. Due 
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Fig. 5: Average value of Pfaffian as a function of Ef, for w = 
20,30,40 and A = 0.64. Inset: comparison of average value of 
Pfaffian for A = 0,M= — Ito 1 — 2 ^ m for same system at 
w = 20. 



to the finite w, even at A = there is a non-zero dis- 
tance between the curves. In the thermodynamic limit, 
the insulator has average Pfaffian ±1, and the metal has 
absolute value of average Pfaffian less than 1. Finite size 
simulations cannot identify the insulating regime with cer- 
tainty, but these simulations suggest that the topological 
insulator, where the index is always —1, is located a short 
distance inside the —2/3 curve and the ordinary insulator 
is located a short distance outside the +2/3 curve. The 
curves with Ef > were obtained by mirroring the curves 
with Ef < 0; while our index is not obviously particle- hole 
symmetric, numeric tests indicate that the average Pfaf- 
fian is in fact close to invariant under this reflection. 

Interestingly, for A ?^ 2, the system has a delocalized 
phase near zero energy, but there is no topological insu- 
lating phase at this value of A and the phase instead sep- 
arates two ordinary insulating phases. Rather, the delo- 
calized phase results from the proximity to the topological 
insulating phase at smaller values of A. 

Relation of Bott Index to Hall Conductance. — 

Finally, we relate the Bott index, in the case without time 
reversal symmetry, to the Hall conductance. We do this by 
relating the Bott index to an invariant [18] which describes 
the response of charge to a magnetic field. Suppose H has 
a spectral gap. Then, we can "spectrally flatten" i7, to 
define a new matrix H' whose eigenvalues are all equal to 
±1, which is still local in the sense that the matrix ele- 
ments of H' are exponentially small between far separated 
sites, and such that the space spanned by negative energy 
states of H^ is the same as that spanned by occupied states 
of H. From now on in this section, when we write H we 
in fact mean the spectrally flattened Hamiltonian H' . 

We add a pseudospin degree of freedom with a Zee- 
mann coupling to a magnetic monopole in the center of 
the sphere, giving the Hamiltonian 
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Coupling strength A 

Fig. 6: For given A, plot of which Ef gives average index equal 
to ±1/3, ±2/3 for k; = 20. 



with the radius of the sphere equal to one so J^^a -^a — 
I. Then, if the system is large, the lattice spacing a is 
small, and we can choose h such that a <^h <^1. Let H 
have M negative eigenvalues. Then H has 2M negative 
eigenvalues since /i <C 1. Since a <^h/\i the pseudospin is 
aligned with the magnetic field, and hops to a nearby site 
within distance a, the angle between the pseudospin and 
the magnetic field on the new site is small. So for a <C 
/i, we can do describe the dynamics semi-classically: the 
particle hops on the lattice, with its pseudospin following 
the magnetic field. In this semi-classical limit, the negative 
eigenvalues are all close to — 1 ± /i. 

The invariant we will consider is the difference between 
the number of negative eigenvalues close to — 1 + /i, which 
we denote M+, and the number of negative eigenvalues 
close to — 1 — /i, which we denote M_, where M+ ± M_ = 
2M. We compute M+ — M_ in two different ways. 

First, we compute it semi-classically. Consider the 
eigenvalues close to — 1 ± /i. If the pseudospin is aligned 
with the field, the Berry phase for transport of the pseu- 
dospin leads to an effective orbital magnetic field, as if 
there were a magnetic monopole inside the sphere. We 
will write P = (1 — H)/2 for the projector onto the neg- 
ative energy states of H. We write P(+) for the projec- 
tor onto the negative energy state of a modified Hamil- 
tonian, i^(±), where H{-\-) is the Hamiltonian H modi- 
fied by an orbital coupling to a magnetic monopole. We 
write P(— ) for the projector onto the negative energy 
states of a Hamiltonian H{—)^ which has the opposite sign 
monopole. So, M^ - M^ = tr(P(±) - P(-)). This quan- 
tity tr(P(±) — P(— )) has been suggested by [18 . 

However, we can also compute M+ — M_ using de- 
generate perturbation theory. Since /i <C 1, we can ig- 
nore matrix elements of h ^^ Xa Ca which couple neg- 
ative and positive energy states of H ^ I and and just 
consider the operator h ^^ Xa cr^ projected onto the 
negative energy space of H ^ L This operator equals 
hj:^PXaP^(Ja = /i-P(PXi,PX2,PX3P). Thus, in 
this semi-classical limit the Bott index equals the differ- 
ence tr(P(±) - P(-))/2. 
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While the invariant of p!8l may be unfamihar, in fact 
physicahy it simply describes the response of charge to 
flux. In a quantum Hall system, the adiabatic addition of 
flux increases the charge by an amount proportional to the 
Hall conductance [19|. It is possible in this way to equate 
the Bott index we consider with the Hall conductance of 
a system on a sphere in the limit of large system size. 

We now turn to the torus. The relation between Hall 
conductance and the Bott invariant as in eq. (H]) will be 
investigated in [7 . The calculation in this section suggests 
that we want a map from the torus to the sphere such 
that the pullback of the connection on CP^ onto T^ gives 
a connection with constant curvature to give a constant 
effective orbital magnetic field on the torus with a single 
monopole. No such map exists but the map we used in 
the log- method, eq. (pTj) . gives constant curvature except 
at the singularity = 0. 

Discussion. — We have presented an approach to cal- 
culating Z2 invariants in time-reversal invariant systems, 
allowing us to study larger system sizes than previously. 

In a numerical study of a model for HgTe, we observed 
a metallic phase around the transition when the spin com- 
ponents are not decoupled, as expected from the existence 
of metallic phases in disordered systems with spin-orbit 
scattering due to anti-localization [24j . 

We use methods of non-commutative topology, not non- 
commutative geometry [20], to define our indices. While 
no precise boundary exists between these fields in the 
mathematics literature, our general approach is to dis- 
card most of the information in the original Hamiltonian, 
only retaining the band-projected position operators. This 
keeps topological properties while discarding "metric in- 
formation" [21^ used in the noncommutative geometry ap- 
proach. By retaining the minimum of information, we 
obtain an efficient, practical computational procedure for 
these indices. The topological invar iance of these indices 
follows from standard linear algebra, although the tighter 
estimates in [9] require recent results in matrix theory. 

Two distinct approaches to investigating topological 
phases are common. One is via index and the other is 
via transport. In this paper we have taken the approach 
of studying an index. In contrast to flux-torus approaches 
to Z2 index, such as Essin-Moore, our approach allows 
the numerical study of much larger systems. For sys- 
tems with coupling between spin components, studying 
the Hall conductance as in [23] does not work, so the 
transport method requires boundaries. In contrast, on 
such systems our method works without the introduction 
of boundaries. Although the transport approach allows 
larger system size, introducing boundaries may increase 
finite size effects. Further, the conductance in a transport 
measurement is not exactly quantized in finite size sys- 
tems, while the index approach returns an integer topo- 
logical invariant even for finite sizes. Another interesting 
ability of the index approach is in studying effects such as 
we have seen near A = 2. There, the delocalized phase 



has a slightly negative average index (roughly —0.3) at 
Ef = 0, implying that the index detects the nearby pres- 
ence of a topological insulator. Ideally, we want access 
to all properties of a topological insulator: bulk conduc- 
tance, conductance of boundary modes, and index. We 
expect that a combination of transport and index tech- 
niques will be needed to obtain a full understanding of 
these systems in the presence of strong disorder. 

Since our approach does not use twisted boundary con- 
ditions, we maintain the time reversal symmetry of the 
problem throughout the calculation. This leads to a philo- 
sophical advantage of our approach: rather than studying 
the response of the system to a perturbation that breaks 
the symmetry (the twisted boundary conditions), we di- 
rectly study an intrinsic property of the given Hamilto- 
nian, in particular allowing us to show obstructions to 
localized Wannier functions for that Hamiltonian. 



MBH thanks M. Freedman, A. Kitaev, C. Nayak and 
X.-L. Qi for useful discussions. 
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